R Notebook

Read in the labelled object Run the correlation comparison between the total scRNA per FACS pop and the AB levels Calculate the proportion of the cell types from the orginal labels UMAP of the cell POPS merged UMAP of the subtypes labelled

library(Seurat)

Attaching SeuratObject
Attaching sp

library(dplyr)

Attaching package: ‘dplyr’

The following objects are masked from ‘package:stats’:

    filter, lag

The following objects are masked from ‘package:base’:

    intersect, setdiff, setequal, union

library('reshape')

Attaching package: ‘reshape’

The following object is masked from ‘package:dplyr’:

    rename

library("ggplot2")

Read in the lableled object single cell sequencing

saveRDS(seu.sc,,paste(pathway, "CombinedLabeledMarkers14102022.RDS"))

Error in saveRDS(seu.sc, , paste(pathway, "CombinedLabeledMarkers14102022.RDS")) : 
  'file' must be non-empty string

NA

Read in the Flow Cytometry data

output_path = "/Users/rhalenathomas/Documents/Data/FlowCytometry/PhenoID/Analysis/PaperFigures/PostSortGatingPops/"
seu.fc <- readRDS(paste(output_path,"june10postSort.seu.RDS", sep = ""))
unique(seu.fc$GatedPop)

[1] Astrocytes RadialGlia Neurons1   Neurons2  
Levels: Astrocytes RadialGlia Neurons1 Neurons2

Marker list

marker.ab <- c("CD24","CD56","CD29","CD15","CD184","CD133","CD71","CD44","GLAST","AQP4","HepaCAM", "CD140a","O4")

marker.genes <-c("CD24","NCAM1","ITGB1","FUT4","CXCR4","PROM1","TFRC","CD44","SLC1A3","AQP4","HEPACAM", "PDGFRA","NKX6-2")

Make the expression matrixes

seu.sc <- ScaleData(object = seu.sc, features = rownames(seu.sc))

Centering and scaling data matrix

  |                                                                                               
  |                                                                                         |   0%
  |                                                                                               
  |===                                                                                      |   3%
  |                                                                                               
  |=====                                                                                    |   6%
  |                                                                                               
  |========                                                                                 |   9%
  |                                                                                               
  |==========                                                                               |  12%
  |                                                                                               
  |=============                                                                            |  15%
  |                                                                                               
  |================                                                                         |  18%
  |                                                                                               
  |==================                                                                       |  21%
  |                                                                                               
  |=====================                                                                    |  24%
  |                                                                                               
  |========================                                                                 |  26%
  |                                                                                               
  |==========================                                                               |  29%
  |                                                                                               
  |=============================                                                            |  32%
  |                                                                                               
  |===============================                                                          |  35%
  |                                                                                               
  |==================================                                                       |  38%
  |                                                                                               
  |=====================================                                                    |  41%
  |                                                                                               
  |=======================================                                                  |  44%
  |                                                                                               
  |==========================================                                               |  47%
  |                                                                                               
  |============================================                                             |  50%
  |                                                                                               
  |===============================================                                          |  53%
  |                                                                                               
  |==================================================                                       |  56%
  |                                                                                               
  |====================================================                                     |  59%
  |                                                                                               
  |=======================================================                                  |  62%
  |                                                                                               
  |==========================================================                               |  65%
  |                                                                                               
  |============================================================                             |  68%
  |                                                                                               
  |===============================================================                          |  71%
  |                                                                                               
  |=================================================================                        |  74%
  |                                                                                               
  |====================================================================                     |  76%
  |                                                                                               
  |=======================================================================                  |  79%
  |                                                                                               
  |=========================================================================                |  82%
  |                                                                                               
  |============================================================================             |  85%
  |                                                                                               
  |===============================================================================          |  88%
  |                                                                                               
  |=================================================================================        |  91%
  |                                                                                               
  |====================================================================================     |  94%
  |                                                                                               
  |======================================================================================   |  97%
  |                                                                                               
  |=========================================================================================| 100%

Calculate correlation between two matrixes

df.cor

           Astrocytes RadialGlia   Neurons1    Neurons2
Astrocytes  0.5553328  0.1522878 -0.4517431 -0.18064326
RadialGlia  0.1305815  0.4451073 -0.3314321  0.05559489
Neurons1   -0.4098105 -0.3016957  0.4534910  0.05175022
Neurons2   -0.4858683 -0.1624758  0.3770386  0.23752196

write.csv(df.cor, paste(output_path,"CorrelationsFCscRNAseq.csv"))

Plot the heatmap Figure 6A

# need to melt the matrix
longData<- melt(df.cor)

Warning in type.convert.default(X[[i]], ...) :
  'as.is' should be specified by the caller; using TRUE
Warning in type.convert.default(X[[i]], ...) :
  'as.is' should be specified by the caller; using TRUE

head(longData)

colnames(longData) <- c("FC","scRNA","R2")

ggplot(longData, aes(x=scRNA,y=FC, fill =R2)) + geom_tile() +
  scale_fill_gradient2(low = "#075AFF",
                       mid = "#FFFFCC",
                       high = "#FF0000") +
   guides(fill = guide_colourbar(label = TRUE,
                                ticks = FALSE)) + theme_bw() +
  coord_fixed()  +scale_x_discrete(expand=c(0,0))+
  scale_y_discrete(expand=c(0,0))+ theme(text = element_text(size=16, colour = "black"),
        axis.text.x = element_text(size = 18, colour = "black", angle = 90), axis.text.y = element_text(size= 18, colour = "black"))  + xlab('scRNAseq') + ylab('Flow Cytometry') 

output_path <- "/Users/rhalenathomas/Documents/Projects_Papers/PhenoID/ForFigures/scRNA/"

pdf(paste(output_path,"FC_scRNA_correlationOct16.pdf"), width = 8, height = 5)

ggplot(longData, aes(x=scRNA,y=FC, fill =R2)) + geom_tile() +
  scale_fill_gradient2(low = "#075AFF",
                       mid = "#FFFFCC",
                       high = "#FF0000", midpoint = 0, limit = c(-0.5,0.5)) +
   guides(fill = guide_colourbar(label = TRUE,
                                ticks = FALSE)) + theme_bw() +
  coord_fixed()  +scale_x_discrete(expand=c(0,0)) +
  scale_y_discrete(expand=c(0,0)) + 
  theme(text = element_text(size=16, colour = "black"),
        axis.text.x = element_text(size = 16, colour = "black", angle = 90), axis.text.y = element_text(size= 16, colour = "black"))  + xlab('scRNAseq') + ylab('Flow Cytometry') 

dev.off()

quartz_off_screen 
                2 

Compare proportions of cells in each FACS pop - Figure 6B

Make a table of cell types

library(stringr)
library(reshape)

# for original cell types 


df.r <- reshape(df, idvar = "Var2", timevar = "Var1", direction = "wide")
dat1 <- df.r
dat1[] <- lapply(dat1[], function(x){
  # Check if the column is numeric
  if (is.numeric(x)){
    return(x/sum(x)*100)
  } else{
    return(x)
  }
})
dat1

write.csv(dat1, paste(output_path,"proportionofCelltypesin4pops.csv"))

### for new main cell types 

df.r <- reshape(df.3, idvar = "Var2", timevar = "Var1", direction = "wide")

Error in reshape(df.3, idvar = "Var2", timevar = "Var1", direction = "wide") : 
  object 'df.3' not found

UMAP with original idents to show overlap (or lack there of) Figure 6C

# Figure 6 C
# UMAP with the 4 orig.idents

sample.order <- c("Astrocytes","RadialGlia","Neurons1","Neurons2")
sample.order <- rev(sample.order)

# colour order to match cell type order
clust.colours <- c("royalblue", "indianred2","palegreen2","springgreen4")
 
Idents(seu.sc) <- 'orig.ident'
      
DimPlot(seu.sc, order = sample.order, cols = clust.colours, shuffle = TRUE, raster=FALSE, pt.size = 0.1, label = FALSE)

output_path <- "/Users/rhalenathomas/Documents/Projects_Papers/PhenoID/ForFigures/scRNA/"

### Figure 6C

pdf(paste(output_path,"UMAPscRNAseqMerge4Oct16.pdf"),width = 7, height = 4)

DimPlot(seu.sc, order = sample.order, cols = clust.colours, shuffle = TRUE, raster=FALSE, pt.size = 0.1, label = FALSE, label.size = 6) +
  theme(legend.text = element_text(size=16), axis.title.y = element_text(size=16), 
        axis.title.x = element_text(size=16), axis.text.y = element_text(size =16),
        axis.text.x = element_text(size =16))
dev.off()

quartz_off_screen 
                2 

## Figure S23 panel 

cell.order <- c("Astrocytes","Endothelial","Glia","NPC","Neurons","Other","Radial Glia")
cell.order <- rev(cell.order)

# colour order to match cell type order
clust.colours <- c("chocolate1","deepskyblue","steelblue4","red2","mediumpurple3","burlywood3", 
                   "pink2")

 
Idents(seu.sc) <- 'Cell_Types'
      
# designated the order of the splits factor
seu.sc$orig.ident <- factor(x = seu.sc$orig.ident, levels = c("Astrocytes","RadialGlia","Neurons1","Neurons2"))


DimPlot(seu.sc, order = cell.order, cols = clust.colours, shuffle = TRUE, raster=FALSE, pt.size = 0.1, label = FALSE, split.by = 'orig.ident', ncol = 2)




pdf(paste(output_path,"UMAP_merge_splitbyorigidentOct16.pdf"),width = 12, height = 7.5)

DimPlot(seu.sc, order = cell.order, cols = clust.colours, shuffle = TRUE, raster=FALSE, pt.size = 0.1, label = FALSE, split.by = 'orig.ident', label.size = 6, ncol = 2) +
  theme(legend.text = element_text(size=16), axis.title.y = element_text(size=16), 
        axis.title.x = element_text(size=16), axis.text.y = element_text(size =16),
        axis.text.x = element_text(size =16))
dev.off()

quartz_off_screen 
                2 

# also plot one not split to see the whole thing


#pdf(paste(output_path,"UMAP_merge_CellTypesOct16.pdf"),width = 6.7, height = 4.1)
DimPlot(seu.sc, order = cell.order, cols = clust.colours, shuffle = TRUE, raster=FALSE, pt.size = 0.1, label = FALSE, label.size = 6) +
  theme(legend.text = element_text(size=16), axis.title.y = element_text(size=16), 
        axis.title.x = element_text(size=16), axis.text.y = element_text(size =16),
        axis.text.x = element_text(size =16))

#dev.off()

Figure 6D - UMAP with subgroups

pdf(paste(output_path,"UMAP_CellSubtypeMarkersLablesOct16.pdf"),width = 12, height = 5)
DimPlot(seu.sc, order = sample.order, cols = clust.colours, shuffle = TRUE, raster=FALSE, pt.size = 0.25, label = TRUE, repel = TRUE)
dev.off()

null device 
          1 

Make a heatmap

Heatmap or Dotplot Save pdf

save the heatmap Change the colours

pdf(paste(output_path,"HeatMapsubtypesOct17.pdf",sep = ""),width = 11, height = 5)
DoHeatmap(seu.sc, features = feature_list, group.by = 'Cell_Subtype_Markers', group.colors = clust.colours, disp.max = 2, disp.min = -1.5,
          angle = 90) + scale_fill_gradientn(colors = c("#154c79", "#eeeee4", "#e28743")) + 
  theme(axis.text.y = element_text(size = 16))

Scale for 'fill' is already present. Adding another scale for 'fill', which will replace the
existing scale.

dev.off()

null device 
          1 

Subtype markers

Idents(seu.sc) <- 'Cell_Subtypes'
neurons <- subset(seu.sc, idents = c("Neurons1","Neurons2","Neurons3","Neurons4",
                                     "Neurons5","Neurons6","Neurons7"))

Idents(neurons) <- 'Cell_Subtype_Markers'

pdf(paste(output_path,"UMAP_NeuronsOct17.pdf",sep = ""),width = 7, height = 4)
DimPlot(neurons) +  theme(text = element_text(size=16, colour = "black"))
dev.off()

null device 
          1 

Idents(seu.sc) <- 'Cell_Subtypes'
neurons.da <- subset(seu.sc, idents = c("DANeurons1","DANeurons2","DANeurons3"))

Idents(neurons.da) <- 'Cell_Subtype_Markers'
pdf(paste(output_path,"UMAP_DANeuronsOct17.pdf",sep = ""),width = 7, height = 4)

DimPlot(neurons.da)+  theme(text = element_text(size=16, colour = "black"))
dev.off()

null device 
          1 

Idents(seu.sc) <- 'Cell_Subtypes'
astro <- subset(seu.sc, idents = c("Astrocytes1","Astrocytes2","Astrocytes3"))
Idents(astro) <- 'Cell_Subtype_Markers'
pdf(paste(output_path,"UMAP_AstroOct17.pdf",sep = ""),width = 7, height = 4)

DimPlot(astro)+  theme(text = element_text(size=16, colour = "black"))
dev.off()

null device 
          1 

Idents(seu.sc) <- 'Cell_Subtypes'
rg <- subset(seu.sc, idents = c("RadialGlia1","RadialGlia2","RadialGlia3","RadialGlia4",
                                "RadialGlia5","RadialGlia6"))
Idents(rg) <- 'Cell_Subtype_Markers'
pdf(paste(output_path,"UMAP_RGOct17.pdf",sep = ""),width = 7, height = 4)

DimPlot(rg)+  theme(text = element_text(size=16, colour = "black"))
dev.off()

null device 
          1 

DotPlots of markers in cell type subsets

n.markers <- c("ASCL1","CP","GRIA2","MGP","SPARCL1","TPH1","TFPI2","CD24")

pdf(paste(output_path,"Dotplot_NeuSubMarkersOct17.pdf",sep = ""),width = 7, height = 4)
DotPlot(neurons, features = n.markers) + 
  theme(text = element_text(size=16, colour = "black"),axis.text.x = element_text(angle = 90))
dev.off()

null device 
          1 

da.markers <- c("RAB3B","TPBG","TTR","TH","SOX6","CALB1","SLC17A6")

pdf(paste(output_path,"Dotplot_DANeuSubMarkersOct17.pdf",sep = ""),width = 7, height = 4)

DotPlot(neurons.da, features = da.markers) + 
  theme(text = element_text(size=16, colour = "black"),axis.text.x = element_text(angle = 90))

Warning: Scaling data with a low number of groups may produce misleading results

dev.off()

null device 
          1 

rg.markers <- c("CYP1B1","NEAT1","PTN","RPL41","TOP2A","VCAN","SOX2","SLIT2","HES1","VIM")

pdf(paste(output_path,"Dotplot_RGMarkersOct17.pdf",sep = ""),width = 7, height = 4)

DotPlot(rg, features = rg.markers) + 
  theme(text = element_text(size=16, colour = "black"),axis.text.x = element_text(angle = 90))
dev.off()

null device 
          1 

# Apoe, Gfap, Aqp4 and Slc1a3
astro.markers <- c("COL3A1","FABP5","HPD","APOE","S100B","IGFBP2","DBI",
                   "PRSS56", "IGTP","LFIT3","LIGP1",
                   "COL1A2")

pdf(paste(output_path,"Dotplot_AstroMarkersOct17.pdf",sep = ""),width = 7, height = 4)

DotPlot(astro, features = astro.markers) + 
  theme(text = element_text(size=16, colour = "black"),axis.text.x = element_text(angle = 90))

Warning in FetchData.Seurat(object = object, vars = features, cells = cells) :
  The following requested variables were not found: IGTP, LFIT3, LIGP1
Warning: Scaling data with a low number of groups may produce misleading results

dev.off()

null device 
          1 

Make a label level of main cell types from the merge data

unique(seu.sc$Cell_Subtype_Markers)

 [1] Neurons-GRIA2     NPC               Neurons-SPARCL1   DANeurons-TPBG    Neurons-CP       
 [6] DANeurons-RAB3B   DANeurons-TTR     RadialGlia-CY1B1  RadialGlia-RPL41  RadialGlia-PTN   
[11] Astrocytes-FABP5  Astrocytes-HPD    RadialGlia-TOP2A  Mix               Neurons-TPH1     
[16] Neurons-TFPI1     Neurons-ASCL1     RadialGlia-NEAT1  RadialGlia-VCAN   Neurons-MGP      
[21] Astrocytes-COL3A1
21 Levels: Astrocytes-HPD Neurons-ASCL1 RadialGlia-RPL41 Neurons-SPARCL1 ... RadialGlia-VCAN

Idents(seu.sc) <- 'Cell_Subtype_Markers'

cluster.ids <- c("Astrocytes","Neurons","RadialGlia","Neurons",
                 "Astrocytes","Neurons","RadialGlia","Neurons",
                 "Astrocytes","Neurons","RadialGlia",
                 "DANeurons","Neurons","Mix","NPC",
                 "Neurons","DANeurons","DANeurons",
                 "RadialGlia","RadialGlia","RadialGlia"
                 )


names(cluster.ids) <- levels(seu.sc)
seu.sc <- RenameIdents(seu.sc, cluster.ids)
seu.sc$Cell_Type2 <- Idents(seu.sc)
Idents(seu.sc) <- 'Cell_Type2'
DimPlot(seu.sc)

DimPlot(seu.sc, group.by = 'Cell_Subtype_Markers')

DimPlot(seu.sc, group.by = 'Cell_Subtypes')

DimPlot(seu.sc, group.by = 'Cell_Type2')

NA
NA

# original clusters proportion of cell types in 
library(reshape2)

pr.celltypes <- as.data.frame(table(seu.sc$orig.ident,seu.sc$Cell_Type2))
pr.celltypes <- reshape(pr.celltypes, idvar = "Var2", timevar = "Var1", direction = "wide")
pr.celltypes

table(seu.sc$Cell_Type2)

Astrocytes    Neurons RadialGlia  DANeurons        Mix        NPC 
      3082       2870       1913        476        228        154 

write.csv(pr.celltypes, paste(output_path,"FreqCellTypes2inFACS.csv"))

Proportion of cell types in original FACS populations

pr

                     1          22         43           64          
Var2                 "Neurons2" "Neurons1" "Astrocytes" "RadialGlia"
Freq.DANeurons-RAB3B "131"      " 76"      "  0"        " 52"       
Freq.DANeurons-TPBG  "10"       "50"       " 7"         "46"        
Freq.DANeurons-TTR   "65"       "20"       " 0"         "19"        

Proportion of cell types tests

library("scProportionTest")

prop_test <- sc_utils(seu.sc)

prop_test <- permutation_test(
    prop_test, cluster_identity = "Cell_Type2",
    sample_1 = "Neurons1", sample_2 = "Neurons2",
    sample_identity = "orig.ident"
)


permutation_plot(prop_test)

prop_test <- permutation_test(
    prop_test, cluster_identity = "Cell_Type2",
    sample_1 = "Neurons1", sample_2 = "RadialGlia",
    sample_identity = "orig.ident"
)
permutation_plot(prop_test)

prop_test <- permutation_test(
    prop_test, cluster_identity = "Cell_Type2",
    sample_1 = "Neurons2", sample_2 = "RadialGlia",
    sample_identity = "orig.ident"
)

*** recursive gc invocation
*** recursive gc invocation
*** recursive gc invocation
*** recursive gc invocation
*** recursive gc invocation
*** recursive gc invocation
*** recursive gc invocation
*** recursive gc invocation
*** recursive gc invocation

permutation_plot(prop_test)

prop_test <- permutation_test(
    prop_test, cluster_identity = "Cell_Type2",
    sample_1 = "Astrocytes", sample_2 = "RadialGlia",
    sample_identity = "orig.ident"
)
permutation_plot(prop_test)

Check the proportions of subtypes of DA neurons

#library("scProportionTest")

prop_test <- sc_utils(neurons.da)

prop_test <- permutation_test(
    prop_test, cluster_identity = "Cell_Subtype_Markers",
    sample_1 = "Neurons1", sample_2 = "Neurons2",
    sample_identity = "orig.ident"
)


permutation_plot(prop_test)

Differential gene expression between Neurons1 and Neurons2 for Other neurons and DA neurons.

down <- rownames(deg.neurons %>% filter(avg_log2FC < -6))
length(down)

[1] 8

In DA neurons

DotPlot(sub.neur, group.by = 'orig.ident', features = up.down) + RotatedAxis()

Warning: Scaling data with a low number of groups may produce misleading results

Look at the GO and other pathway analysis for the DGE

For DA neurons

pdf(paste(output_path,"GOresultesNeurons1vs2.pdf"), width = 12, height = 6)
plotEnrich(t.GObio.da1, showTerms = 10, numChar = 40, y = "Count", orderBy = "Combined.Score") +
  theme(text = element_text(size=16, colour = "black"))
plotEnrich(t.GObio.da1, showTerms = 10, numChar = 40, y = "Count", orderBy = "Overlap")+
  theme(text = element_text(size=16, colour = "black"))
plotEnrich(t.GObio.da2, showTerms = 10, numChar = 40, y = "Count", orderBy = "Combined.Score")+
  theme(text = element_text(size=16, colour = "black"))
plotEnrich(t.GObio.da2, showTerms = 10, numChar = 40, y = "Count", orderBy = "Overlap")+
  theme(text = element_text(size=16, colour = "black"))
dev.off()

null device 
          1 

Dot plots of some up regulated genes

DotPlot(all.neurons, group.by = 'orig.ident', features = regulate.genes) + RotatedAxis()

Warning: Scaling data with a low number of groups may produce misleading results

pdf(paste(output_path,"DotPlotDEG_N1vsN2.pdf"))

DotPlot(all.neurons, group.by = 'orig.ident', features = regulate.genes) + RotatedAxis()

Warning: Scaling data with a low number of groups may produce misleading results

dev.off()

quartz_off_screen 
                2 

BiocManager::install("clusterProfiler")

'getOption("repos")' replaces Bioconductor standard repositories, see '?repositories' for details

replacement repositories:
    CRAN: https://cran.rstudio.com/

Bioconductor version 3.15 (BiocManager 1.30.18), R 4.2.1 (2022-06-23)
Installing package(s) 'clusterProfiler'
also installing the dependencies ‘blob’, ‘plogr’, ‘Biostrings’, ‘fastmatch’, ‘ggfun’, ‘ggplotify’, ‘tidytree’, ‘treeio’, ‘DBI’, ‘RSQLite’, ‘KEGGREST’, ‘DO.db’, ‘fgsea’, ‘aplot’, ‘scatterpie’, ‘shadowtext’, ‘ggtree’, ‘AnnotationDbi’, ‘downloader’, ‘DOSE’, ‘enrichplot’, ‘GO.db’, ‘GOSemSim’, ‘qvalue’, ‘yulab.utils’

trying URL 'https://cran.rstudio.com/bin/macosx/contrib/4.2/blob_1.2.3.tgz'
Content type 'application/x-gzip' length 46035 bytes (44 KB)
==================================================
downloaded 44 KB

trying URL 'https://cran.rstudio.com/bin/macosx/contrib/4.2/plogr_0.2.0.tgz'
Content type 'application/x-gzip' length 13193 bytes (12 KB)
==================================================
downloaded 12 KB

trying URL 'https://bioconductor.org/packages/3.15/bioc/bin/macosx/contrib/4.2/Biostrings_2.64.1.tgz'
Content type 'application/x-gzip' length 14356457 bytes (13.7 MB)
==================================================
downloaded 13.7 MB

trying URL 'https://cran.rstudio.com/bin/macosx/contrib/4.2/fastmatch_1.1-3.tgz'
Content type 'application/x-gzip' length 49267 bytes (48 KB)
==================================================
downloaded 48 KB

trying URL 'https://cran.rstudio.com/bin/macosx/contrib/4.2/ggfun_0.0.7.tgz'
Content type 'application/x-gzip' length 193540 bytes (189 KB)
==================================================
downloaded 189 KB

trying URL 'https://cran.rstudio.com/bin/macosx/contrib/4.2/ggplotify_0.1.0.tgz'
Content type 'application/x-gzip' length 137271 bytes (134 KB)
==================================================
downloaded 134 KB

trying URL 'https://cran.rstudio.com/bin/macosx/contrib/4.2/tidytree_0.4.1.tgz'
Content type 'application/x-gzip' length 249644 bytes (243 KB)
==================================================
downloaded 243 KB

trying URL 'https://bioconductor.org/packages/3.15/bioc/bin/macosx/contrib/4.2/treeio_1.20.2.tgz'
Content type 'application/x-gzip' length 916609 bytes (895 KB)
==================================================
downloaded 895 KB

trying URL 'https://cran.rstudio.com/bin/macosx/contrib/4.2/DBI_1.1.3.tgz'
Content type 'application/x-gzip' length 745805 bytes (728 KB)
==================================================
downloaded 728 KB

trying URL 'https://cran.rstudio.com/bin/macosx/contrib/4.2/RSQLite_2.2.18.tgz'
Content type 'application/x-gzip' length 4509316 bytes (4.3 MB)
==================================================
downloaded 4.3 MB

trying URL 'https://bioconductor.org/packages/3.15/bioc/bin/macosx/contrib/4.2/KEGGREST_1.36.3.tgz'
Content type 'application/x-gzip' length 185051 bytes (180 KB)
==================================================
downloaded 180 KB

trying URL 'https://bioconductor.org/packages/3.15/bioc/bin/macosx/contrib/4.2/fgsea_1.22.0.tgz'
Content type 'application/x-gzip' length 1489172 bytes (1.4 MB)
==================================================
downloaded 1.4 MB

trying URL 'https://cran.rstudio.com/bin/macosx/contrib/4.2/aplot_0.1.8.tgz'
Content type 'application/x-gzip' length 56478 bytes (55 KB)
==================================================
downloaded 55 KB

trying URL 'https://cran.rstudio.com/bin/macosx/contrib/4.2/scatterpie_0.1.8.tgz'
Content type 'application/x-gzip' length 404495 bytes (395 KB)
==================================================
downloaded 395 KB

trying URL 'https://cran.rstudio.com/bin/macosx/contrib/4.2/shadowtext_0.1.2.tgz'
Content type 'application/x-gzip' length 226478 bytes (221 KB)
==================================================
downloaded 221 KB

trying URL 'https://bioconductor.org/packages/3.15/bioc/bin/macosx/contrib/4.2/ggtree_3.4.4.tgz'
Content type 'application/x-gzip' length 916106 bytes (894 KB)
==================================================
downloaded 894 KB

trying URL 'https://bioconductor.org/packages/3.15/bioc/bin/macosx/contrib/4.2/AnnotationDbi_1.58.0.tgz'
Content type 'application/x-gzip' length 5175190 bytes (4.9 MB)
==================================================
downloaded 4.9 MB

trying URL 'https://cran.rstudio.com/bin/macosx/contrib/4.2/downloader_0.4.tgz'
Content type 'application/x-gzip' length 21727 bytes (21 KB)
==================================================
downloaded 21 KB

trying URL 'https://bioconductor.org/packages/3.15/bioc/bin/macosx/contrib/4.2/DOSE_3.22.1.tgz'
Content type 'application/x-gzip' length 6704314 bytes (6.4 MB)
==================================================
downloaded 6.4 MB

trying URL 'https://bioconductor.org/packages/3.15/bioc/bin/macosx/contrib/4.2/enrichplot_1.16.2.tgz'
Content type 'application/x-gzip' length 283459 bytes (276 KB)
==================================================
downloaded 276 KB

trying URL 'https://bioconductor.org/packages/3.15/bioc/bin/macosx/contrib/4.2/GOSemSim_2.22.0.tgz'
Content type 'application/x-gzip' length 920638 bytes (899 KB)
==================================================
downloaded 899 KB

trying URL 'https://bioconductor.org/packages/3.15/bioc/bin/macosx/contrib/4.2/qvalue_2.28.0.tgz'
Content type 'application/x-gzip' length 2801783 bytes (2.7 MB)
==================================================
downloaded 2.7 MB

trying URL 'https://cran.rstudio.com/bin/macosx/contrib/4.2/yulab.utils_0.0.5.tgz'
Content type 'application/x-gzip' length 33694 bytes (32 KB)
==================================================
downloaded 32 KB

trying URL 'https://bioconductor.org/packages/3.15/bioc/bin/macosx/contrib/4.2/clusterProfiler_4.4.4.tgz'
Content type 'application/x-gzip' length 821345 bytes (802 KB)
==================================================
downloaded 802 KB

The downloaded binary packages are in
    /var/folders/k4/khtkczkd5tn732ftjpwgtr240000gn/T//RtmpOFXQrS/downloaded_packages

installing the source packages ‘DO.db’, ‘GO.db’

trying URL 'https://bioconductor.org/packages/3.15/data/annotation/src/contrib/DO.db_2.9.tar.gz'
Content type 'application/x-gzip' length 1769978 bytes (1.7 MB)
==================================================
downloaded 1.7 MB

trying URL 'https://bioconductor.org/packages/3.15/data/annotation/src/contrib/GO.db_3.15.0.tar.gz'
Content type 'application/x-gzip' length 29908485 bytes (28.5 MB)
==================================================
downloaded 28.5 MB

* installing *source* package ‘DO.db’ ...
** using staged installation
** R
** inst
** byte-compile and prepare package for lazy loading
** help
*** installing help indices
** building package indices
** testing if installed package can be loaded from temporary location
** testing if installed package can be loaded from final location
** testing if installed package keeps a record of temporary installation path
* DONE (DO.db)
* installing *source* package ‘GO.db’ ...
** using staged installation
** R
** inst
** byte-compile and prepare package for lazy loading
** help
*** installing help indices
** building package indices
** testing if installed package can be loaded from temporary location
** testing if installed package can be loaded from final location
** testing if installed package keeps a record of temporary installation path
* DONE (GO.db)

The downloaded source packages are in
    ‘/private/var/folders/k4/khtkczkd5tn732ftjpwgtr240000gn/T/RtmpOFXQrS/downloaded_packages’
Old packages: 'BiocParallel', 'brew', 'cli', 'clue', 'commonmark', 'cpp11', 'crayon', 'curl',
  'data.table', 'DelayedMatrixStats', 'devtools', 'digest', 'dotCall64', 'evaluate', 'foreign',
  'gert', 'ggforce', 'ggraph', 'ggridges', 'graphlayouts', 'igraph', 'irlba', 'isoband',
  'jsonlite', 'lifecycle', 'limma', 'mnormt', 'nlme', 'nnet', 'openssl', 'pracma', 'purrr',
  'ragg', 'rainbow', 'RcppArmadillo', 'RCurl', 'readr', 'rlang', 'rmarkdown', 'sctransform',
  'Seurat', 'SeuratObject', 'SingleCellExperiment', 'sys', 'testthat', 'tidyselect', 'tinytex',
  'vctrs', 'vroom', 'XML', 'yaml'

Update all/some/none? [a/s/n]: 

n
BiocManager::install("pathview")

'getOption("repos")' replaces Bioconductor standard repositories, see '?repositories' for details

replacement repositories:
    CRAN: https://cran.rstudio.com/

Bioconductor version 3.15 (BiocManager 1.30.18), R 4.2.1 (2022-06-23)
Installing package(s) 'pathview'
also installing the dependencies ‘KEGGgraph’, ‘org.Hs.eg.db’

trying URL 'https://bioconductor.org/packages/3.15/bioc/bin/macosx/contrib/4.2/KEGGgraph_1.56.0.tgz'
Content type 'application/x-gzip' length 1665739 bytes (1.6 MB)
==================================================
downloaded 1.6 MB

trying URL 'https://bioconductor.org/packages/3.15/bioc/bin/macosx/contrib/4.2/pathview_1.36.1.tgz'
Content type 'application/x-gzip' length 2700402 bytes (2.6 MB)
==================================================
downloaded 2.6 MB

The downloaded binary packages are in
    /var/folders/k4/khtkczkd5tn732ftjpwgtr240000gn/T//RtmpOFXQrS/downloaded_packages

installing the source package ‘org.Hs.eg.db’

trying URL 'https://bioconductor.org/packages/3.15/data/annotation/src/contrib/org.Hs.eg.db_3.15.0.tar.gz'
Content type 'application/x-gzip' length 83788492 bytes (79.9 MB)
==================================================
downloaded 79.9 MB

* installing *source* package ‘org.Hs.eg.db’ ...
** using staged installation
** R
** inst
** byte-compile and prepare package for lazy loading
** help
*** installing help indices
** building package indices
** testing if installed package can be loaded from temporary location
** testing if installed package can be loaded from final location
** testing if installed package keeps a record of temporary installation path
* DONE (org.Hs.eg.db)

The downloaded source packages are in
    ‘/private/var/folders/k4/khtkczkd5tn732ftjpwgtr240000gn/T/RtmpOFXQrS/downloaded_packages’
Old packages: 'BiocParallel', 'brew', 'cli', 'clue', 'commonmark', 'cpp11', 'crayon', 'curl',
  'data.table', 'DelayedMatrixStats', 'devtools', 'digest', 'dotCall64', 'evaluate', 'foreign',
  'gert', 'ggforce', 'ggraph', 'ggridges', 'graphlayouts', 'igraph', 'irlba', 'isoband',
  'jsonlite', 'lifecycle', 'limma', 'mnormt', 'nlme', 'nnet', 'openssl', 'pracma', 'purrr',
  'ragg', 'rainbow', 'RcppArmadillo', 'RCurl', 'readr', 'rlang', 'rmarkdown', 'sctransform',
  'Seurat', 'SeuratObject', 'SingleCellExperiment', 'sys', 'testthat', 'tidyselect', 'tinytex',
  'vctrs', 'vroom', 'XML', 'yaml'

Update all/some/none? [a/s/n]: 

n
BiocManager::install("enrichplot")

'getOption("repos")' replaces Bioconductor standard repositories, see '?repositories' for details

replacement repositories:
    CRAN: https://cran.rstudio.com/

Bioconductor version 3.15 (BiocManager 1.30.18), R 4.2.1 (2022-06-23)
Warning: package(s) not installed when version(s) same as current; use `force = TRUE` to re-install:
  'enrichplot'
Old packages: 'BiocParallel', 'brew', 'cli', 'clue', 'commonmark', 'cpp11', 'crayon', 'curl',
  'data.table', 'DelayedMatrixStats', 'devtools', 'digest', 'dotCall64', 'evaluate', 'foreign',
  'gert', 'ggforce', 'ggraph', 'ggridges', 'graphlayouts', 'igraph', 'irlba', 'isoband',
  'jsonlite', 'lifecycle', 'limma', 'mnormt', 'nlme', 'nnet', 'openssl', 'pracma', 'purrr',
  'ragg', 'rainbow', 'RcppArmadillo', 'RCurl', 'readr', 'rlang', 'rmarkdown', 'sctransform',
  'Seurat', 'SeuratObject', 'SingleCellExperiment', 'sys', 'testthat', 'tidyselect', 'tinytex',
  'vctrs', 'vroom', 'XML', 'yaml'

Update all/some/none? [a/s/n]: 

n
library(clusterProfiler)

clusterProfiler v4.4.4  For help: https://yulab-smu.top/biomedical-knowledge-mining-book/

If you use clusterProfiler in published research, please cite:
T Wu, E Hu, S Xu, M Chen, P Guo, Z Dai, T Feng, L Zhou, W Tang, L Zhan, X Fu, S Liu, X Bo, and G Yu. clusterProfiler 4.0: A universal enrichment tool for interpreting omics data. The Innovation. 2021, 2(3):100141

Attaching package: ‘clusterProfiler’

The following object is masked from ‘package:reshape’:

    rename

The following object is masked from ‘package:stats’:

    filter

library(enrichplot)
# we use ggplot2 to add x axis labels (ex: ridgeplot)
library(ggplot2)

Try with cluster profiler

emapplot(gse, showCategory = 10)

Error in has_pairsim(x) : 
  Term similarity matrix not available. Please use pairwise_termsim function to deal with the results of enrichment analysis.

Predict cell types again to make table 12 Organoids in house data

seu.r <- AIW120
colnames(seu.r@meta.data)

 [1] "orig.ident"                         "nCount_RNA"                        
 [3] "nFeature_RNA"                       "percent.mt"                        
 [5] "predicted.id"                       "prediction.score.RGa"              
 [7] "prediction.score.Astrocytes.2"      "prediction.score.Neurons.e"        
 [9] "prediction.score.Oligodendrocytes"  "prediction.score.Neurons.DA"       
[11] "prediction.score.Neural.Precursors" "prediction.score.Neurons.i"        
[13] "prediction.score.Neurons"           "prediction.score.Other"            
[15] "prediction.score.Epithelial"        "prediction.score.Astrocytes.1"     
[17] "prediction.score.RGd1"              "prediction.score.RGd2"             
[19] "prediction.score.max"               "pred.165days"                      
[21] "RNA_snn_res.0"                      "RNA_snn_res.0.05"                  
[23] "RNA_snn_res.0.25"                   "RNA_snn_res.0.5"                   
[25] "RNA_snn_res.0.8"                    "seurat_clusters"                   
[27] "RNA_snn_res.0.9"                    "RNA_snn_res.1"                     
[29] "RNA_snn_res.1.2"                    "res08names"                        
[31] "res08names2"                        "res08names.gene2"                  
[33] "main.genotype"                      "res08names.gene3"                  
[35] "res08names3"                       

anchors <- FindTransferAnchors(reference = seu.r, query = seu.q, dims = 1:25)

Performing PCA on the provided reference using 2000 features as input.
Projecting cell embeddings
Finding neighborhoods
Finding anchors
    Found 4595 anchors
Filtering anchors
    Retained 2761 anchors
as(<ngCMatrix>, "dgCMatrix") is deprecated since Matrix 1.5-0; do as(., "dMatrix") instead

print("getting predictions")

[1] "getting predictions"

predictions <- TransferData(anchorset = anchors, refdata = seu.r$res08names2, k.weight = 50)

Finding integration vectors
Finding integration vector weights
0%   10   20   30   40   50   60   70   80   90   100%
[----|----|----|----|----|----|----|----|----|----|
**************************************************|
Predicting cell labels

seu.q <- AddMetaData(seu.q, predictions$predicted.id, col.name = "prediction")
seu.q <- AddMetaData(seu.q, predictions$prediction.score.max, col.name = "prediction.score.max")


seu.q$AIW120.pred <- ifelse(seu.q$prediction.score.max > 0.8, seu.q$prediction, "none")
DimPlot(seu.q, group.by = 'AIW120.pred')

More predictions Developing brain: Cortex Forebrain

Adult brain: whole brain with subtypes Adult brain midbrain and striatum Adult brain midbrain main cell types Adult brain DA or Astro subtypes

Tables of top predictions

t.lables <- as.data.frame(table(seu.q$Cell_Subtype_Markers,seu.q$devcortex))
t.lables$Freq <- as.double(t.lables$Freq)

top.prediction <-as.data.frame(t.lables  %>% group_by(Var1)  %>% top_n(1, Freq))
top.prediction